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Abstract 

When matter falls past the horizon of a large black hole, the expectation from string 
theory is that the configuration thermalizes and the information in the probe is rather 
quickly scrambled away. The traditional view of a classical unique spacetime near a 
black hole horizon conflicts with this picture. The question then arises as to what 
spacetime does the probe actually see as it crosses a horizon, and how does the back- 
ground geometry imprint its signature onto the thermal properties of the probe. In this 
work, we explore these questions through an extensive series of numerical simulations 
of DO branes. We determine that the DO branes quickly settle into an incompressible 
symmetric state - thermalized within a few oscillations through a process driven en- 
tirely by internal non-linear dynamics. Surprisingly, thermal background fluctuations 
play no role in this mechanism. Signatures of the background fields in this thermal 
state arise either through fluxes, i.e. black hole hair; or if the probe expands to the size 
of the horizon - which we see evidence of. We determine simple scaling relations for the 
DO branes' equilibrium size, time to thermalize, lifetime, and temperature in terms of 
their number, initial energy, and the background fields. Our results are consistent with 
the conjecture that black holes are the fastest scramblers as seen by Matrix theory. 
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1 Introduction and results 



The correspondence between strongly coupled gauge theories and quantum gravity p]-[l] 
suggests that gravitation and spacetime may be viewed as emergent structures - emergent 
from within the intricate non-linear interactions of non-gravitational quantum field theories. 
One aspect of this correspondence that remains particularly challenging to understand has 
to do with a thermodynamics theme that seems to underly gravitational dynamics [21 [5] . 

A representative illustration of this problem arises in the following process. A stringy - 
possibly Planckian - probe falls into a black hole. The expectation from the gravity-gauge 
theory correspondence is that the probe gets thermalized as it flies through the spacetime 
region near the horizon of the black hole [6] . What underlies this mechanism of thermalization 
on the gravity side? More precisely, what is the spacetime that the probe actually sees as 
it crosses the horizon - if not the smooth, traditional, no-hair region that ends with a 
pathological singularity. Furthermore, there have been suggestions that this thermalization 
process is an unusual one, characterized by a fast scrambling of the information in the 
probe [7]. The situation also ties in with the black hole information paradox that attempts 
to account for such scrambled information [91 Dili EH]. In a slightly different language, we 
want to find out why does the traditional picture of the spacetime near a black hole horizon 
fail to fully capture horizon physics - even when the curvature scales at the horizon are very 
small. 

In this work, we investigate these questions using the strongly coupled + 1 dimensional 
gauge theory that describes DO brane dynamics p21 021 H] • We want to quantify the ther- 
malization process on the gauge theory side using numerical simulations; and identify the 
role played by background fields that the DO branes are immersed in. If a spacetime - be 
it that of a traditional black hole or one of the many fuzzball geometries [HI [T3] - is to 
thermalize a probe, in what way does this spacetime fix the thermodynamic attributes of 
this thermalization process? For example, if a probe is to get scrambled into a configuration 
of temperature determined by the size of a black hole horizon, we want to find out how 
the geometrical information about the size of the horizon gets eventually encoded into the 
temperature attribute of the scrambled probe. Through this, we can start addressing the 
difficult question as to how one determines the emergent spacetime that an in-falling probe 
actually experiences - assuming the traditional no-hair geometry somehow falters in the 
vicinity of a horizon. 

The numerical simulation of + 1 dimensional U(N) Super Yang-Mills theory at strong 
coupling is a difficult one - even without the inclusion of supersymmetry — since it is par- 
ticularly computationally intensive (TBJ HZl [TH1 UH]- Larger values of N provide for numerical 
stability at the expense of time. And a thermodynamic treatment requires a large ensemble 
of simulations. We overcome these problems by a series of efficient physical and technical 
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tricks that we develop; and through the use of recent technological advances in parallel pro- 
cessing. The result is a framework of DO brane dynamics exploration that can be done in 
real-time, pocking and tweaking the parameters as the dynamics evolves so as to develop 
physical intuition about this very rich system. 

We track the time evolution of a probe consisting of N DO branes in various background 
field configurations - in their center of mass frame. The DO brane coordinates are represented 
by N x N hermitian traceless matrices nine in total, one for every spatial direction. The 
eigenvalues of these matrices may be interpreted as the positions of the individual DO branes. 
We define the extent of the probe in a subspace of the nine space directions as 

pl=Tr(xV) (1) 

where /i is one representative index within the subspace of interest (no sum over /i). The 
actual size of the probe is defined as [20] 

for jj, summed over 1 ... 9. In all cases, we start the DO branes in a static two-sphere 
configuration in three of the nine space directions 

x i = s a 1 , [a 1 , a j ] = ie ijk a k For i, j, k = 1,2,3 (3) 

where the a*'s are the Pauli matrices in an N x iV representation, and where s is a tunable 
parameter that fixes the initial size. This initial size is then given by 

R = s V2VN 2 - 1 ~ s V2N (4) 

for iV ^> 1. Tuning s also corresponds to adjusting the initial energy of the setup: smaller s 
corresponds to smaller initial energy. We refer to the six directions transverse to the initial 
spherical configuration as the transverse directions; and we refer to the remaining three as 
probe directions. Matrices in the three probe directions are denoted by x % with % = 1,2,3; 
while those in the transverse directions are denoted by y a with a = 1 ... 6, initially set to 
zero. The momenta matrices in the probe directions are denoted by p\ while those in the 
transverse directions are labeled q a . We also add small random gaussian fluctuations to the 
initial positions of the DO branes - off the spherical shape of the probe. This is so as to avoid 
starting on top of a saddle point in the potential, ending up exploring non-generic regions 
of the phase space. This initial randomization of the positions can be viewed as accounting 
for the effect of initial quantum or thermal fluctuations [21]. 
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Figure 1: The evolution of DO branes in zero but fluctuating background fields. In this case, 
we have 15 DO branes whose initial energy is much larger than the energy in the fluctuating 
background fields. The horizontal axis is time r in all three graphs. The smaller two graphs 
are plots of a selection of eigenvalues of the DO branes: in one of the three probe directions 
on the left, and in one of the six transverse directions on the right. In the middle, the graph 
shows the extent of the probe: in red and blue for two probe directions, and in green for a 
transverse direction. All variables are dimensionless, as defined in the main text. 
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In this work, we study over 500 simulations, all qualified with a numerical error of the 
order of 10%. Figure [I] is a stereotypical time evolution when the initial conditions are 
perturbed slightly. Without initial perturbations, we find oscillatory evolution with no ther- 
malization. Throughout, lengths are measured in units of £, and energies and temperatures 
in units of l/£, where 

t= (2vr) 2/3 ^I » Is , (5) 
9.' 

with the string coupling g s <C 1 and the string length given by l s . The spherical configuration 
oscillates to a smaller extent in the probe directions; at the same time, its initial zero extent 
in the transverse directions expands to form a uniform fuzzy ball in all nine space directions. 
In a short time - within a few oscillations - the eigenvalue spectrum becomes thermal. The 
effective Yang-Mills coupling is mostly large - i. e. the thermalization phenomenon is a strong 
coupling effectj^j There are various interesting attributes of this process that we discover as 
we vary the background parameters. We explore three qualitatively different scenarios. 
First, we consider zero but thermally fluctuating background fields. This corresponds to 
exploring the effects of a dilute gas of massless supergravity fields onto the probe. In this 
analysis, we include the effects of thermal back-reaction - a central but delicate mechanism 
to the process of thermalization in traditional thermodynamic systems. We consider only a 
fluctuating background metric and D2 brane flux for simplicity. Next, we consider non-zero 
background fields with negligible fluctuations and study the effects of non-zero gravity and 
D2 brane flux on DO probe thermalization. Finally, we model the problem of the probe 
falling into a large Schwarzschild black hole and determine the effect of the in-fall on the 
thermalization of the probe. We summarize the results in the upcoming three sections. 

1.1 Results: zero backgrounds 

We start with zero background fields - with or without thermal fluctuations of the super- 
gravity fields. Even when the background fluctuations are zero, we introduce a level of 
randomization into the problem by perturbing the initial spherical configuration of the DO 
branes at around 10% — 50% level. We can vary N, the number of DO branes in the probe, 
and the initial size of the probe sphere - through the parameter s which also tunes the ini- 
tial energy. Figure [2] shows the data we collect from a simulation for s = 1 and iV = 15 DO 
branes. In this case, we have no background field fluctuations. These six plots are generated 
for every simulation, and we use this example to outline the typical analysis we perform 

3 By this we mean that the non-linear terms in the Hamiltonian evolution equations play a central role. 
The dynamics is classical chaotic only if these terms are initially large, and we find that, if the initial energy 
is tuned such that the non-linearities are weak, there is no thermalization. 
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for any simulation. The top three plots in the figure show the eigenvalues and extent of 
some of the position matrices. The lower three plots attempt to quantify the thermalization 
phenomenon. The two leftmost graphs on the bottom are histograms of the eigenvalues of 
momenta matrices, computed for the last half of the simulation timeline where the probe has 
settled into a potentially thermal configuration. A semicircle pattern is indicative of ther- 
malization, with the radius of the semicircle giving us a measure of the probe's temperature 
t p = r\y) |2T]. The rightmost plot on the bottom is a graph of the autocorrelation function 



as a function of time r; and the O operator is given by O = Tr (X 1 + iY 1 ) . This measures 
the autocorrelation of a gauge invariant operator over time: for a thermalizing state, the 
function should decay exponential over time. We see from the example depicted in the 
figure that the probe does indeed thermalize. From these graphs, we extract the timescale 
of thermalization r t h and the probe temperature t p as shown. Eventually, our simulation 
invariably breaks down. We track this by looking at the constraint equation 



which needs to be satisfied for a consistent dynamical system. This break-down is seen to 
occur often (but not always) because the size of the probe suddenly explodes. It signals the 
accumulation of significant numerical errors by the evolution algorithm. For zero background 
fields, there is a flat direction in the system when all nine matrices become mutually com- 
muting. Then the configuration can expand unbounded with no cost of energy - the larger 
in size the better from entropic considerations. It is likely that the numerical errors that we 
see in the simulations are due to this instability some of the time. Hence, we associate a 
lifetime with the probe, denoted by tj, the time it takes for the probe to de-stabilize and fly 
apart, measured from the start of the simulation as shown in the Figure. We are uncertain 
whether this quantity T\ is a physical one, or a pathology of numerical simulations. Finally, 
from a blow-up of the size plot, we can read off the equilibrium size of the probe R eq . 

Hence, each simulation gives us four outputs: r t h, u, t p , and R eq . We want to determine 
the dependence of these quantities on three input parameters: the initial size of the spherical 
probe tuned by s - equivalent to keeping track of the energy; the number of DO branes N; and 
when we turn on thermal background field fluctuations, the temperature t of the background 
gas. 

Figure [3] shows a compilation of the results from a series of simulations which fix all 
parameters and scan over s. Figure [4] shows a compilation of the results from simulations 
which fix all parameters and scan over N instead, the number of DO branes. For all of these 
simulations, fluctuations of the background fields are turned off. If the background fields are 



Re<C t (0)C(r)) 



(6) 



[x\ Pi ]+[y a ,q a } = 



(7) 
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Figure 2: A sample simulation, along with a depiction of the output parameters: r t h for the 
time to thermalize, T; for lifetime, t p = r\ h for probe temperature, and R eq = a/9 x p x »,/\fN 
for final probe size. For this simulation, the input parameters were: N = 15, s — 1, with 
no background fluctuations. In the probe extent and distribution graphs, red represent the 
probe direction x 1 , blue is for the probe direction x 2 , and green is for the transverse direction 

y l - 
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Figure 3: Results from a collection of simulations showing the output parameters r t h (ther- 
malization time), T[ (lifetime), R eg (final equilibrium size), and t p = r\ h (probe temperature) 
as a function of the input parameter s. For all simulations, N = 15, and there are no 
background fluctuations. All slopes are within the indicated nearest rational values that is 
allowed by the 10% error of the simulations. 
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Figure 4: Results from a collection of simulations showing the output parameters r t h (ther- 
malization time), n (lifetime), R eq (final equilibrium size), and t p = rf ft (probe temperature) 
as a function of the input parameter N, the number of DO branes. For all simulations except 
the ones used to gauge the lifetime 7j, s = 0.03, with no background fluctuations. For the 
graph of the lifetime dependence on N, another set of simulations were used due to technical 
limitations with storing large matrices for long simulations. For the latter case, s = 0.5 
with much shorter simulation times, with no background fluctuations, and larger statistical 
errors. All slopes are within the indicated nearest rational values that is allowed by the 10% 
error of the simulations. 
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made to thermally fluctuate with temperature t, there is no signature of this temperature 
in the probe's thermalization dynamics as depicted in Figure |5j We have scanned over 
several orders of magnitude of background temperatures t within the regime of validity of 
our formalism, and considered the delicate effect of thermal back-reaction - with no change 
in the conclusion. Putting all these results together, we find the following scaling relation 
for the equilibrium size R eq 



For the thermalization time, we find 



R eq oc s N 



IN 

Tth oc - oc — . (9) 



And for the temperature of the probe t p , we find 

i? 4 



r 



th -t p * , 4 iV 2 oc . (10) 
As for the lifetime of the probe, we find 

11 

T l <X Ar q/9 °C = . (11) 

sN 3 / 2 R eq y/N 

Irrespective of the details, the DO brane probe seems to eventually disintegrate in a time T\. 
The thermalization time r t h is of particular interest given the proposal of [7J |8] suggesting 
that black holes are highly efficient scramblers - with their scrambling timescale proportional 
to the logarithm of the number of degrees of freedom. We will comment on the implications 
of our results with regards to this proposal in the Conclusion section at the end. 

If the probe is too small and the background gas temperature t is too large, the probe's 
thermalization is disrupted and the evolution follows deterministic oscillations. We cannot 
determine whether this is a physical effect or a numerical pathology since these simulations lie 
at the edge of the regime of validity of our formalism. We also consider 'on-shell background 
fluctuations': that is, background field fluctuations that satisfy the Laplace equation at the 
center of mass of the probe. This corresponds to a scenario where there is no matter sourcing 
the background fields at the probe's location. Once again, thermalization is disrupted and 
we are unable to determine whether the effect is physical or numerical. Either way, these 
scenarios are uninteresting since they do not lead to thermalization. 



All four of these relations (|8j) to ( 11 ) are insensitive to the temperature of the background 
t\ A thermal bath of massless supergravity fields seems to play no role in the thermalization 
of the probe. 
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Figure 5: The effect of background ffuctuations. For all graphs, we have t = 1, N = 15, 
s — 5, and g s = 10~ 4 . We see that background fluctuations do not effect the thermal 
properties of the probe, except for its lifetime. The size of the fluctuations is determined by 
the temperature t, as well as the string coupling g s , as described in detail in the main text. 
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The conclusive theme of the analysis can be summarized as follows: the thermalization of 
a probe of DO branes is an internal process arising from the non-linear DO brane interactions 
and strong Yang-Mills coupling. A quartic term in the DO brane matrix coordinates in the 
Hamiltonian of the system underlies this thermalization phenomenon. The result depends 
on the number of DO branes and their equilibrium size. Fluctuations in the background 
fields do not play a relevant role in this process, a fact that is rather counter-intuitive as 
compared to more traditional statistical mechanical systems. There is however a role played 
by the fluctuating fields in destabilizing the configuration at high enough temperatures. This 
regime is however at the edge of the regime of validity of our analysis. 

The size of the probe may be sensitive to large non-zero background fields, as opposed 
to small fluctuating ones. We will investigate this case next. 

1.2 Results: quasi-static non-zero backgrounds 

We next consider non-zero background fields with negligible thermal fluctuations. We turn 
on the effect of the background metric through tidal gravitational forces acting on the probe 
in its center of mass frame; and we also turn on background D2 brane flux. These fields 
may be arising from a large number of background D-branes that our probe is inserted into. 
If the background is to represent a black hole, we may expect D2 brane flux from previous 
investigations of Matrix black holes (221 [23J El] . 

We track the dimensionless curvature scales in the three and six dimensional subspaces 
through the two parameters mi and respectively. As Figure [6] shows, we find no effect of 
mi and vri2 on lifetime, thermalization time, temperature, or the final size of the probe. That 
is unless the length scale associated with the m's becomes small enough to compete with 
the size of the probe. This is akin to the gravitational Gregory-Laflamme phenomenon, and 
syncs well with the gravitational analogue of this scenario presented in |25j. At this point, 
we get a ringing behavior as shown in Figure [7j Thermalization is disrupted and we do 
see a complex correlation between the background's and the probe's thermal characteristics. 
Unfortunately, this regime tests the bounds of validity of our simulation and cannot be 
reliably quantified yet. 

We track the dimensionless D2 brane flux through a parameter labeled n. We do find a 
correlation between n and the thermal properties of the probe as shown in Figure [8] We are 
able to identify the following scaling laws 

R eq ocn , r t h oc n~ 1/2 ; (12) 

r th = t P oc n 2 , 77 oc n~ . (13) 

These results indicate a sensitivity of the thermal properties of the probe to black hole hair, 
if present. Larger D2 flux results in larger thermal equilibrium size of the probe, which 
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mi = 2 , m 2 = -2 m x = 0.2 , rrt 2 = -0.2 m x = 0.02 , m 2 = -0.02 m 1 = 0.002 , m 2 = -0.002 




Figure 6: The effect of scanning over the gravitational tidal force parameters m\ and 777.2. 
A selection of simulation are shown, but many more were analyzed with similar qualitative 
conclusions. For all simulations in this figure, s = 0.5, iV = 15, with no background fluctu- 
ations. At around \m,i\ = |m 2 | = 2, we see the start of the probe size competing with the 
length scale associated with the background. 
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Probe extent p 2 xll X 1 eigenvalues P',P 2 distributions 

Figure 7: When the probe size becomes bigger that the length scale characterizing the 
background metric, we see a ringing effect and the break-down of thermalization. The 
validity of our simulations is also breaking down in this region of the parameter space. For 
this simulation, we have m\ = m,2 =, n = —0.05, s = 0.5, iV = 15, with no background 
fluctuations. 

correlates with the proximity of a larger background black hole. Larger D2 brane flux seems 
to also hasten the thermalization. The probe's equilibrium temperature seems however 
to increase with the background flux, a rather counter-intuitive scaling relation from the 
perspective of the naive Matrix black hole model [221 123| El]. 

1.3 Results: the in-fall problem 

In the third and last scenario, we drop the probe DO branes from rest near the vicinity of the 
Schwarzschild black hole. We start at an initial radial distance ro from the horizon such 
that ro/r/j = 10; and we track the evolution as a function of local Fermi normal coordinate 
time, as the probe crosses the horizon and is sucked into the black hole singularity. We 
consider a large black hole to assure that our simulation are reliable near the horizon. In 
units of l/£, we choose the temperature of the hole as t = 0.004. In units of £, the flight 
time to the singularity is Tfi ~ 336, while the time to reach the horizon is Th, or ~ 325. And 
Thor is well within the regime of validity of our simulation. Figure [9] summarizes the results 
for a collection of iV = 15 DO branes. We find that the probe thermalizes due to internal 
dynamics as seen earlier, when it is away from the horizon. We find surprisingly that this 
thermalization is insensitive to the black hole temperature! The thermalization time and the 
equilibrium probe temperature and size do not depend on the black hole parameters. As the 
probe approaches the horizon, we do however see the probe explode in size. Our formalism 
however breaks down before the probe can reach the size of the horizon. The trend however 
is consistent with a probe expanding violently in size as it nears the horizon, as proposed 
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Figure 8: Results from a collection of simulations showing the output parameters (ther- 
malization time), n (lifetime), R eq (final equilibrium size), and t p = r^ h (probe temperature) 
as a function of the input parameter n, the D2 brane flux . For all simulations, N = 15, and 
there are no background fluctuations. 
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Figure 9: The probe's evolution as it falls into a Schwarzschild black hole. The black hole 
temperature is t bh = 0.004, corresponding to tidal force parameter mi = 6 x 10~ 7 as defined 
in the main text. The other parameters of the simulation are N = 15 for number of DO 
branes, g s = 10~ 6 for the string coupling, and r /rh = 10. In the top two graph, we plot the 
probe's size and the probe's radial coordinate distance from the black hole as a function of 
time. In the lower left figure, we see the eigenvalue evolution and the thermalization prior 
to reaching the horizon. Finally, in the lower right graph, we plot the e parameter, defined 
later on in the main text: this parameter tracks the adiabatic regime of validity of the Fermi 
normal coordinates. If e ~ 1, the simulation results cannot be trusted. As we can see, the 
parameter remains well within the reliable range as the probe crosses the horizon, but not 
much beyond. 
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Ill 

We leave the conclusions of these results to the Conclusion and Outlook Section at the 
end, Section [6j Before that however, we present the details progressively. In Section [2j we 
develop the problem and lay out the parameter space we want to explore. In Section [3j the 
Langevin problem is formulated for our Matrix theory system, including a delicate thermal 
back-reaction mechanism. In Section |4| the regime of validity of the three cases we study is 
elaborated in detail. And in Section [5j some of the technical numerical considerations are 
presented. 

2 The setup 

Consider N DO branes in arbitrary supergravity background fields. The non-Abelian action 
is given in terms of the coordinates of the DO branes - represented by N x iV U (N) matrices 
$^ with /I = 1 • • • 9 



S = -\ f dtSTr{e- 4 '{-{P[E m + E^{Q- 1 -5Y u E v0 ])) ll2 {Y)etQ) 

Qs^s J 



1/2 



g s i s 



(14) 



where A = 2irl 2 s , = 5% + iA[$ M , & p ]E pu , and E pv = G^ v + B pv - all given in the notation 
used in [20] ■ This action is derived using T-duality symmetry. As an expansion in powers of 
A, it is known to agree with direct open string computations to order A 3 |26j. Our goal is to 
study dynamics of DO branes in a certain restricted class of backgrounds to order A 2 . Our 
assumptions follows: 

1. We focus on a non-trivial background metric G pv and 3- form RR potential a 
constant dilaton = g s , and set all other supergravity fields to zero. We want to 
analyze the evolution of a probe consisting of N DO branes falling through background 
fields that may represent a Schwarzschild black hole. Hence, we need to at least keep 
track of the background metric seen by the probe. In the Matrix black hole picture, 
there are hints of a D2 brane structure stretched at the would-be black hole horizon [23] . 
Hence, we also consider a non-zero flux from C^. For simplicity, we eliminate all other 
fields. 

2. We consider an expansion of this action to order A 2 . This results in restrictions on the 
size of the probe, its speed of evolution, and the background metric. We will list these 
restrictions later. 
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3. At order A 2 , it is easy to check that the U(l) sector of the coordinates $ M - corre- 
sponding to the center of mass of the probe - decouples from the SU (N) sector. We 
focus on the relative dynamics of the probe DO branes and drop the U(l) center of 
mass dynamics. Put differently, we will track the evolution of the probe DO branes in 
their center of mass frame. 

4. In the center of mass frame of the probe, we adopt Fermi normal coordinates. This 
means that the metric at the origin is Minkowski, its first derivatives are zero, and the 
second derivatives of the metric - the tidal forces - provide the leading gravitational 
effects. We also assume that the probe does not back-react onto the background, except 
thermally - as we will explain in detail later. 

Under these conditions, we are then left with the Hamiltonian [25] 



H = ± — f- Tr <^ + M^<W + iN aM [<£>", $ 7 - - $"] 2 \ (15) 

where the $ M 's are N x N hermitian traceless matrices representing the non-commutative 
coordinates of N probe DO branes. Having chosen the static gauge for the gauge field on the 
world-line, we also need to supplement the equations of motion with the constraint 

^,$ M ]=0. (16) 

The physical position coordinates of the DO branes in nine space directions are given by 

X» = 2nl 2 ^ . (17) 



The background metric appears through^] 



= —Gtt,^ ] (18) 



while the D2 brane flux appears through 



No* = \ C Z,,] ■ (19) 



4 M A1J/ also depends on DO brane flux = dC^ x \ as shown in [27]. Hence, our analysis throughout 
can be thought of including the effects of DO brane flux in the background through the consideration of the 
M^ v . The existence of non-zero C*- 1 - 1 however also introduces a magnetic, velocity dependent term in the 
Hamiltonian that we have dropped. This implies that our analysis does include the leading effect of RR 
flux for small DO brane speeds. 
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See [27J for a more general and detailed derivation of these relations. 

It is convenient to use dimensionless variables instead - labeled as m^, n a p 7 , and r 
- defined as follows 

W=\P , M, u = ^f , N aPl = ^f , t = rt, (20) 

where £ has dimension of length. This leads to the Hamiltonian 

H = Tr |i(^) 2 + m^P + in aM [<p a , <P P ] <F ~ \ [P, P?\ (21) 

where derivatives are now with respect to dimensionless time r. We will choose i so that 
the scale of the energy is given by l/£ 

&S£.i i (22) 

9.1. < 4 t { 1 

which implies 



£ = (2tt) 2 / 3 ^ • (23) 
Hence, all lengths/times are measured in units of £, and energies in l/£. We then can write^] 

= iTr|^+m^0^ + m Q/37 [0 a ,0T0^-^r,0l 2 | (25) 

where e is dimensionless energy. Denoting by [s] the numerical scale in the dimensionless 
matrices M 

Scale of matrix entries in 0^ ~ [s] =>• Probe size ex [s] £ . (26) 

Since we will need g s <C 1, this means that £ ^> l s . Hence, we are looking at large probes 
and small energies compared to the string scale. 



5 An alternative scaling used in the literature is given by the choice of I such that 

(2^) 2 1 _ 1 
g s l s £^ Qs^s 

i.e the DO brane mass, which implies £ — (27r) 1 / 2 / s , the string scale. 



(24) 
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We define the effective dimensionless coupling in the theory as 

, _g* _ (2tt)-V7 3 _ 1 

This is the parameter that tunes the strength of the quartic and cubic terms in the Hamil- 
tonian. For strong coupling effects, we would expect g 2 ^ ~ 1. For large iV 3> 1, we may 
encountered an effective coupling given by N/e 3 . The simulations presented in this work 
with clear thermalization effects correspond to g\^ > 1. As a rule of thumb, when the di- 
mensionless scale of the size of the probe goes below s < 0.3, we get into the strong coupling 
regime. Note that classically this strong coupling regime is meant as a statement about the 
relative importance of the non-linear (cubic and quartic) terms in the Hamiltonian. 



2.1 Regime of validity 



The expansion of the original Dirac-Born-Infeld (DBI) action ( 14 ) is valid under the condi- 
tions 

[ S ]^/ 3 «l , [s] 2 [m] gt /3 « 1 , [% s 2/3 «l, (28) 

where the square brackets signify 'numerical scale of, [m] denotes the typical scale of the 
background fields m^, and [n] that of n a p y . These three statements follow from the con- 
vergence of the expansion of the DBI. In addition, if the probe size is small relative to the 
length scales in the background fields 

[s} 2 [m] < 1 , [s] 2 [n] 2 < 1 , (29) 

we then expect that the background fields will evolve slowly compared to the evolution of the 
probe degrees of freedom. We do not impose this last restriction when studying background 
fluctuations, but we will need it for the remaining analysis. 

As for the background fields, first the string coupling must be small so that the leading 
supergravity regime is valid 

9s < 1 ■ (30) 
We also need weak curvature scales in the background to protect from excited string states 

[M] << C 2 ^ M <7 S 2/3 << 1 , (31) 
and similarly weak fluxes for the background field n 

M g 1 / 3 « i • (32) 
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The energy of a configuration will have the following scaling structure 

[e]^([m}[s] 2 + [n][sf + [s} 4 )xN\ (33) 

where we have written 

[Tr[(j) 2 }] ~ iV 7 (34) 

with 7 being typically a number between one and three and N being the number of DO 
branes. For example, 7 = 3 for the highly ordered spherical DO brane configuration satisfying 
the SU(2) algebra. As the system evolves through the classical equations of motion, we 
expect the configuration will be attracted toward a state where all three terms in the energy 



expression (33) are of the same order 



[m] [s] 2 ~ [n] [sf ~ [s] 4 [s] 2 ~ [m] , [s] ~ [n] . (35) 

In summary, we have a parameter space consisting of the string coupling g s , the back- 
ground fields rri/xu and n a ^, the DO initial size scale s, and the number N of DO branes. When 
we consider background thermal fluctuations, the scales of and n Q/ g 7 get set through the 
input temperature t of the heat bath as described later. The following conditions delineate 
this parameter space: 



Small string coupling 
Weak supergravity fields 



9s < 1 • (36) 



[m]^ 3 «l , [n]^/ 3 «l. (37) 



Valid DBI expansion 

[ S ]^/ 3 «l , [ S ] 2 [m]^ 3 «l , [% s 2/3 «l. (38) 

Small probe (for non-fluctuating backgrounds) 

[s] 2 [m] < 1 , [s] 2 [n] 2 < 1 . (39) 



To simplify the problem further, we will look at background fields with a lot of symmetry. 
Even with background fluctuations, we expect that the most symmetric scenario is the closest 
to equilibrium. Hence, restricting from the outset to symmetric setups should not miss the 
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final equilibrium states we are interested in. We divide the nine dimensional space in two 
subspaces 

ft-tx* , <p a ->y a (40) 

with = 1 • • • 3, a,b, . . . = 4 • • • 9. The momentum canonical to x % is labeled p l , and 

that to y a is labeled q a . We then consider backgrounds with SO (3) x SO (6) symmetry. This 
means we have two background fields in ra^ v 

fi^u 

m l _ LU ={ mi [i — v — 1,2,3 (41) 
m 2 fM — v — 3, • • • , 9 

which we call m,\ and m^, and one field in n a ^ 

_ f ne al37 a, 0,7 = 1,2,3 , 
Ha ^ ~\0 Otherwise 1 > 

which we simply call n. In part of the discussion, we will freeze the six dimensional subspace 
and consider probe dynamics in 3 + 1 dimensions only. In this setting, we split m x into m 
and mi as follows: 

m /i = v = 1 

"V = { m i At = ^ = 2,3 (43) 
oo Otherwise (dynamics frozen by hand) 

for SO (2) cylindrical symmetry in three dimensions. This is so as to capture the scenario 
of a probe falling into a black hole: the 1 direction is the in-falling radial direction, while 2 
and 3 are transverse to it. 

For the BMN model [28] . we note the special background field values in our notation 

m x = m2 = i , n = \ ■ ( 44 ) 
Our system however encompasses a much larger class of dynamics which we will next explore. 



3 Fluctuating backgrounds 

Part of our program is to determine the role of fluctuating background fields on the process 
of thermalizing probe DO branes. These fluctuations can be quantum mechanical in nature 
or due to a thermal bath. Either way, the problem involves a separation of timescales: 
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one slow time scale characterizing the probe evolution, and another fast one characterizing 
the background fluctuations. A key ingredient in determining the correct evolution of the 
probe is the effect of back-reaction: the probe changes the background fluctuations, which 
in turn modifies the probe's dynamics. This is a delicate mechanism that can be studied for 
quantum fluctuations, as well as thermal fluctuations. We focus in this work on fluctuations 
of thermal nature and will come to the quantum fluctuation case in a future work. To model 
the effect of thermal back-reaction for our system of DO branes, we start with a simpler toy 
example that shares similarities to our case, and expand from there. 

3.1 A toy example 

The process of thermalization in general relies on an interesting back-reaction mechanism. 
In short, a system undergoing thermalization by being in contact with a thermal reservoir is 
associated with two time scales, one much shorter than the other one (see, for example [29]). 
The shorter timescale comes from the fluctuations in the thermal bath. The longer timescale 
is emergent from thermalization and characterizes the evolution of the system towards ther- 
malization. The fluctuations of the thermal bath kick the system by exchanging energy with 
it in such a way that the system's energy is not conserved over the longer timescales. This 
shifts the energy balance in the thermal bath by a relatively small amount, which however 
in turn can affect the bath's fluctuations significantly enough to react back on the system 
differently. This delicate mechanism is at the heart of thermalization in numerous physical 
systems and underlies for example Brownian motion dynamics. To illustrate it concretely, 
and in particular as it is relevant to our problem of DO brane dynamics, consider a system 
consisting of a single particle of mass m subject to a spring force with energy 

E = -mu 2 + -kr 2 . (45) 

In this setup however, the spring constant k is not constant and is fluctuating ergodically 
around an average value k 

k = (k) th + 5 k = k + 5 k (46) 

where the th subscript on the averaging indicates a thermodynamic ensemble average in the 
heat bath system. The thermal fluctuations Sok are due to the thermal reservoir having a 
fixed temperature T and occur on a short timescale 5ot ~ 1/T. Note that the thermal average 
is then, according to the ergodic theorem of thermodynamics, the same as the time average 
over timescales much larger than Sot. These are attributes of thermodynamic equilibrium in 
heat baths. 

We want to understand how the system evolves under the influence of these fluctuations 
in k and hence how the system may thermalize. Denote the longer timescale associated with 
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the evolution of the system by St, and St 3> S t by assumption. We define a time averaging 
scheme for observables based on this longer timescale as 

/t+st 
dt'f(t') (47) 

for any function of time f{t). The equations of motion will then naively become 

m(a) s = -{kr) s ~ -{k) th (r)s - -k{r)s , (48) 

where the ergodic theorem and St 3> 5 t was used. This would be the end of the story if 
back-reaction is not taken into account: we say the reservoir is so big that it is unaffected 
by our system, and the spring constant is independently maintained at its average value 
k. Our system then oscillates harmonically with average spring constant k. This treatment 
would miss the interesting effect of thermalization that, depending on the circumstances, 
may affect the dynamics significantly. Instead, we need to keep track of the effect of the 
evolution of the system on the heat bath which then back reacts on the system. 

The fluctuations Sok occur over the short timescale Sot and have a prescribed distribution 
determined by the nature of the thermal ensemble. For example, for high enough temper- 
atures and low enough densities, we may take the fluctuations in k to obey the classical 
Maxwellian distribution at equilibrium 

(<W 2 

Prob e ,(5 A;) = e -^= (49) 

where the standard deviation a k can be related to the details of the heat reservoir such as its 
temperature. The fluctuations in Sok must not be correlated over time scales greater than 
Sot. We would expect the correlation function to decay fast as in 

(S k(t) Sok{t + St)) th ^e' St /^ . (50) 

Focus on an instant in time t', with t < t' < t + St, the interval over which the position of 
the particle changes appreciably. To leading order, the thermal average would be 

(S k(t')) th ~ (51) 

if we were to ignore back-reaction effects. Over the time interval S't — t' — t, the system's 
parameters have changed by 

v ->• v + S'v =>• E ->■ E + S'E . (52) 
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In particular, 5'E is not zero since the energy exchange of our system with the heat bath 
would change the system's energy. In the process, the heat bath's energy changes also by an 
amount —5'E. Hence, we may write immediately 

Prob(t', 5 k) _ _ ps , E 

Prob(t, 5 k) ' K ' 



where (3 = 1/T and Proh(t, 5 k) = Prob eq (5 k) is the equilibrium probability given by (49). 
This implies that 

(5 Q k(t'))th ^ (5 k(t)e^ s ' E ) th ~ (5 k(t)(l - 0S'E)) th = -P(5 k(t)S'E) th (54) 
for small [55' E. We can easily compute 5'E as 

5'E = f dt"5 k(t") r(t") ■ 6 -4^- = f dt"5 k(t") r(t") ■ v(t") (55) 
Jt ot" J t 

i.e. energy is conserved if 5ok(t") = 0. Putting things together, one then gets 
(5 k(t')} th ~ -p f dt"r(t")-v(t") (5 k(t)5 k(t")) th 



~ -P r{t) ■ v(t) (5 k(0) 2 ) th 5 t = -(3 r{t) ■ v{t) a\ 5 t (56) 

where in the last steps we used the fact that (5ok(t)5 k(t")) th is significant only over the 
range Sot, and that r(t) and v(t) change little during this interval. Note that this expression 
is now non-zero and leads to a velocity dependent correction to the force law, i.e. leads to 
a dissipation mechanism. To see this, we look back at the equation of motion (48), and now 
write 

nt+St 

m(a) s ~ -ik) th (r) s - / dt"5 k(t")r(t") ~ -k(r) 5 - (5 k(t)) th (r) 5 



= -k(r) s + p<Tl6 t(r)s(r)s-{v)8 = -k(r) s + a-(v) s (57) 

with the tensor 

a = p 4 5 t(r)s(r) 5 . (58) 

This is the Langevin equation for the system, and a is the dissipation tensor. The effect of a 
is crucial in thermalizing the spring system. The timescale for thermalization is then easily 
identified as 

St ~ — > 5 t . (59) 
m 
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This hierarchy between timescales, fast fluctuations driving a thermalization mechanism on 
a slower timescale, is the subtlety underlying the approach to equilibrium for most systems. 

There is a quicker way to reach at this conclusion that will be useful for us in simulating 
Matrix black hole dynamics later on. We may have directly jumped through the sequence 
of approximation by proposing that the probability distribution of the fluctuations S k is 
modified by the presence of the ensemble to 

Prob(t, S k) = Piob eq {S k) x e~^ E = PTob eq (S k)e-^° k r Sr (60) 

where Prob eq (Sok) is the equilibrium probability distribution at temperature 1/(3. For the 
classical limit of fluctuations which this problem is valid for, this leads to 

Prob(t, 5 k) oc e^e-V 5 ° k r - 5r . (61) 

Or equivalently the average value of the fluctuations to leading order shifts from zero to 

(5ok) th -> ~/3r ■ 5ra 2 k (62) 

with the effect of the thermal bath back-reacting onto the spring system included. This 
means that the fluctuating force on the spring is 

- (S k) th r ->■ (3 r • Sr a\ r = (3 a\ S t r (r • v) (63) 

since the fluctuations occur over a timescale Sot. This is again what we obtained earlier 



in (57). The difference is that we have folded the back-reaction effect as a leading order 



contribution by changing the distribution of fluctuations according to the energetics of the 



heat bath. The advantage of using this method, as given by (60), will become apparent 
when we attempt to implement this effect in our numerical simulation. Otherwise, it would 
be hopelessly computationally time intensive to simulate this dynamics with large enough 
Matrix theory matrices. 

3.2 Thermal bath and Matrix theory 



For the case of our Matrix theory Hamiltonian given by (25), we will track fluctuations in 



and n ai a 7 . Adopting the prescription in (60), we then write 



Prob(t, <5 7»v) = Prob eg (5 "v) x exp(-/3 Sm^ 5 [Tr^^]) (64) 

for the fluctuation distribution of m^. Here, 5q [Tr ^^> u \ measures then the change in the 
DO brane matrix coordinates under an evolutionary time step of Sot. Note that (3 is the 
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inverse of the dimensionless temperature; i.e, we are writing energies and temperatures in 
units l/£. For n a( g 7 , we then have 

Prob(t, 5 n a ^) = Prob eg (£ n Q( 3 7 ) x exp(-/3 5n aM S [STr i 7 ]) . (65) 

The question remains as to how to determine the equilibrium distributions Probe^o"^) and 



Prob e(? (5on Qj g 7 ) of the background fields. From (18) and (19), we know that these fluctuations 



correspond to massless excitations in the supergravity fields of thermal origin. We then need 
to look at the supergravity background to determine the energy content of a disturbance 
Sm ^ and Sn a ^. We look at the energy content in the supergravity part of the combined DO 
branes-background system because we expect that the energy is overwhelmingly stored in 
the background fields. That is, we assume the probe is small enough that it can be treated 
as a probe of the background field configuration. 

In our conventions, the supergravity action is given by 



(66) 



where the RR field strengths are F^ = dC« and = dC®, while the NSNS field 
strength is = dB^ 2 \ The dilaton background is taken constant as discussed earlier. 
We are interested in general scaling behavior of the fluctuations with respect to tempera- 
ture and string coupling. The eventual goal is to simulate the fluctuations and study the 
response of the system as a function of temperature and string coupling - as opposed to 
actually according physical numerical importance to, for example, the temperature value. 
We are hence looking for extracting the qualitative features of the distribution functions 
Piob eq (Som^ u ) and Prob e g(5o^a/3 7 ) from supergravity. We can read off the energy content of 
a fluctuation by looking at the time-time component of the energy momentum tensor in the 
Einstein frame, which is still given by 

T a = \ (R ab - (l/2)G ab R) (67) 

K 

since the Einstein tensor is unchanged under frame change involving a constant dilaton, and 
where 2 k 2 = (27r) 7 /^. From the equations of motion, this gives for example 

Ttt= (27M2^ + "' (68) 

for the term coming from the electric field part of the 2-form field strength. Writing the 
field strengths in Taylor expansion about the center of mass of the DO branes, we would get 
contributions of the form 

Fu ^ C t ,i 1 + C m 1 x 3 ->■ C m 1 x 3 (69) 
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where the linear term is irrelevant except for the center of mass motion of the configuration: 
this is in the U(l) part of the DO matrices, which entirely decouples from the SU(N) as can 
be checked from the form of our action. The term M M „ in the DBI includes contributions 
from Cfl as shown in [27]. Hence, we have F ti ~ 5M x L, where this relation is to be read 



as relating the scale of F ti to typical scale of our background fields bM^ and the length scale 
L over which the background fields vary. 

For the energy content of graviton fluctuations, we use Wald's prescription [30] for iden- 
tifying energy with the graviton field 

Tu = ^ (rII } - (1/2)G% } RW) (70) 

where the {2} superscript indicate expansion of the corresponding quantity to second order 
in metric perturbations h^ u about the local flat metric 

R { £ - (i/2)^ } ^ 2 > ~ VA M + ■ • • ■ (Ti) 

Once again, Taylor expanding about the center of mass of the DO branes, we have the 
identification of scales 

h tt ,i - Vclo + htttilo^ > h tt ,ij\ x j ~ 5M x L (72) 



using (18). One then repeats this scale analysis for every supergravity field. Putting things 
together, one then gets for any RR field fluctuations a relation to the associated typically 
energy scale SE 

5E ~ l\- D L D+2 5M 2 ~ l\- D L D 5N 2 ; (73) 
While for NSNS fields, we get an additional factor of g 2 

6E~\ l]- D L D+2 5M 2 ~ \ l]- D L D 5N 2 . (74) 
9 s 9 S 



5M and 5N denote the scales of the background fields in (18) and (19); L is the spatial size 
of the excitation; and D is the number of non-compact dimensions, taken as nine or three. 
The rest of the directions of space are assumed to be compact of size set by the string scale 
l s , most simply through toroidal compactification. 

First, note that for fixed energy scale SE, the size of fluctuations 6M and 6N is dominated 



by the RR fields since we also need j, < 1. Hence, we need to consider only (73). The 
background fields all consist of massless degrees of freedom. In a thermal state, we then 
expect a black body spectrum. The modes 5M and 5N are in general functions of frequency 
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it), but they are peaked around the peak of the black body spectrum u p k ~ T, with a width 
in frequency of the order of the temperature as well 5u ~ T. The length scale L is the 
characteristic scale over which our classical background fields vary. In general, an estimate 
of the number n 7 of massless particles in a box of size L is given by 



"- -T D L D ^L~nll D l (75) 



where 1/T is the thermal wavelength at temperature T. We also know that SE ~ n^uj p k ~ 
n 7 T. One may think that perhaps we need n 7 ^> 1 to have a reliable description of the 
background spectrum through classical fields. This is not necessarily true. For example, 
one can easily describe a coherent state with less than unity occupancy using classical fields. 
However, in our case, we will want n 7 ^> 1: as we shall see, we will adopt a distribution of 
fluctuations that does not capture black body physics for low occupancy numbers; this is done 
for computational efficiency In terms of dimensionless parameters, including temperature 
t written in units of l/£, we then end up with estimates of the fluctuation sizes for the 
background fields as 

8 m ~ sf" 1 )/ 6 t^ D+3 ^ 2 n~ l ' D , 5n ~ gf-^ 6 t^ D+ ^ 2 . (76) 

For example, the size of the fluctuations of the electric field in a gas a photons in a black 
body thermal configuration in three space dimensions would scale as T 3 . This is one power of 
temperature more than the standard 5E 2 ~ T 4 scaling (where E is the electric field) because 
our probe is sensitive to the gradient of the electric field about its center of mass. Hence 
we get an additional power of the thermal wavelength which goes as 1/T. In the same way, 
the gravitational effects on the probe DO branes are tidal in nature and do not involve the 
net gravitational force on its center of mass. The end result is summarized with the simple 



scaling relations shown in (76). 



We take the equilibrium distribution function for the thermal background of massless 
fields as Gaussian 



exp 

Prob eg (5m) = 



Sm 2 
'2al 



exp 



Sn 2 

Prob e9 (H = j==^-, (T7) 



for each mode 5m^ v and 5n Q/ g 7 , with the standard deviations <r^ and a 2 given by (76) 

a m = gf- 1 ^ 6 t^ D+ ^ 2 n- 1 ' 2D , a n = gf~^ t^' 2 . (7* 
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We now use equality signs by absorbing any order unity numerical coefficients in our estimates 
into the definition of our temperature parameter t and coupling g s . There is however one 
additional point we need to be particularly careful about. Given that this is a fluctuation 
near the peak frequency u p k ~ T, we need to apply the fluctuation at every time interval 
1/T in our simulation. Hence, we can use this distribution as long as we perturb the system 
at a frequency near the peak of the black body spectrum. 

Note that this distribution does not capture all quantum aspects of the black body 
distribution. To see this, the energy in, say, a background massless field 5m{u) at frequency 
cu scales as 5m{u) 2 ; and this is related to the average number n 7 (u;) of massless particles of 
type 5m by 

5m 2 (u))du = a(co)n.y(uj)duj = a(co) r— -dco (79) 

where oe(u) is some proportionality factor scaling as uo D ~ l x u>. Hence, the square of the 
classical field measures the number of massless particles at a given frequency. The hallmark 
of the black body spectrum is that the fluctuations in the number of photons are not small 
even for n 7 (w) ^> 1, unlike the ideal gas distribution. This can be expressed as 



An 2 1 



7 



— 2 — 



1 (80) 



for the black body spectrum. Contrast this with classical Maxwellian distribution fluctua- 
tions which obey 

AN 2 1 , , 



N 2 N 

which goes to zero for large N 3> 1. Hence, to capture the correct quantum statistics, we 
need to have 



(dm 2 - 5m 2 } = a 2 An 2 = a% 7 (1 + n 7 ) = 5m 2 (l + 5m 2 } (82) 



for every u, and in particular for the peak u p k ~ T. The gaussian distribution (77) gives 
instead 

5m 2 - 5m 2 } = 2 (jrrA . (83) 

Thus, our fluctuations do remain important for large n 7 3> 1 - as needed from a black body 
spectrum - but not for n 7 of order one or less. In short, treating the background fields 
classically, we need to assume a large condensate of massless particles within a size given 
classical field profiles. Our distribution however will fail to capture the statistics for low 
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particle occupancy. This regime also corresponds to fluctuations that have a parametrically 
small effect on the probe DO branes. Hence, for the purposes of tracing the thermalization 
evolution of our probe, the normal distribution we use is very much adequate. This trick is 
needed to improve the computational efficiency of our simulations. It will translate into a 
lower bound on temperature that we need to restrict our simulations to, as we shall see. 



4 Parameter space 



In this section, we combine the general validity bounds determined in section |2.1| with the 
three particular cases we want to focus on. For each case, we want to determine the regime 
of the available input parameters g s , t, m, n, s, and t we can trust. The three cases are: 

• Case I: DO probe in background fields that average to zero, but otherwise fluctuate 
thermally. This scenario corresponds to the probe embedded in a dilute gas of massless 
supergravity fluctuations. 

• Case II: DO probe in background fields that average to non-zero values, with negligible 
thermal fluctuations. This scenario will attempt to explore correlations of quasi-static 
background fields with the thermal properties of the probe. 

• Case III: DO probe falling into a Schwarzschild black hole with no thermal fluctuations. 
This scenario aims at determining whether the approach of the probe to a black hole 
horizon has a special role in thermalizing the probe. In particular, we would want to 
see whether the probe is prevented to reach the center of the black hole through fast 
thermalization. 

We next proceed in analyzing the 'interesting' regimes in the parameter space for each 
of these cases: a little qualitative analytical control over the dynamics will help in zeroing 
onto the highlights of the numerical simulations without exploring uninteresting deserts of 
parameter values or physically invalid regimes. 

4.1 Case I: Zero fields with thermal fluctuations 

We take the background fields m,\ and m 2 (henceforth collectively referred to as m for 



simplicity), and n averaging to zero with Gaussian distribution given by (77). The scale of 



the background fields is set by the standard deviations (78) and/or the average back-reaction 



shifts obtained from (64) and (65) 



{5 back m}~°^{s} 2 N^ , [5 back n}^^[ s fN^ . (84) 
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Given that a m and a n are determined by the temperature of the background gas, we are thus 
trading the scale of [m] and [n] for the scale of the temperature [t]. Hence, our parameter 
space in this case consists of g s , s, t, and N. 

The interesting mechanism to explore in this scenario has to do with the effect of thermal 
back-reaction. This becomes important when 

cr m < 8 th m , a n < 6 th n . (85) 

We would like to identify the relevant values of g s , s, t, and N for which (1) thermal back- 



reaction is important and (2) the validity conditions outlined in Section 2.1 are satisfied. 



For these purposes, we can focus on fluctuations in m or n since they lead to the same 
conclusions. We know that 

a m ~ g s b t * , (86) 

and 



Equation (85) then becomes 



Sbackm ~ N 1 s 2 g s 3 t 1 



1 > g s > N D-1 S D-H 3 D-1 



The condition of weak background fields from (37) translates to 



where we combined the statement with the weak string coupling condition g s 1 from (|36j) . 

m < g-'^n ^t<^g s 3(V+2) S ~^N~^T2 . (89) 



While the DBI expansion gives us the conditions from (38) 

[ S ]^ 3 «l , [5] 2 [m]^ 3 «l. (90) 
Finally, since our statistical distribution does not capture the case of fluctuations of a small 



number of quanta, we get from ( 75 ) and n 7 ^> 1 

ts > 1 



(91) 



These four curves are depicted in Figure [10| We then restrict our input parameters g s , s, t, 
and iV such that we are within the depicted shaded region in this Figure. 



31 



Figure 10: The parameter space delineating the regime where the simulations can be trusted 
for Case I. For this plot, we have taken g s = 10~ 4 , 7 = 2, N = 15, D = 9, and n 7 = 10 2 . 
There are three thick lines that bound the relevant shaded region in between: a bound arising 
from requiring small background curvatures (blue curve), another from convergence of the 
DBI action expansion (red curve), and another arising from statistical considerations (black 
curve, with n 7 ~ 100). The thin orange line bounds the region where thermal back-reaction 
starts becoming important. All other bounds are weaker that the ones depicted in the input 
parameter range we consider. 
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4.2 Case II: Quasi-static backgrounds 



In the second scenario, we consider classical non-zero background field condensates, on top of 
which we may add small thermal fluctuations. That is, the equilibrium distributions become 



Prob eg (5ra) 
Prob eq (Sn) 



exp 


(6m— m c ) 2 


[ \ 




V 27l(7 m 


exp 


(5n-n c ) 2 





(92) 



where m c and n c are classical background field values fixed by hand. For simplicity, we also 
make sure that the thermal fluctuations do not overwhelm the average values of the back- 
ground fields. In this case, our parameter space is given by g s , s, m, n, and N. Temperature 
t will then be inconsequential. 

To assure that we are within the regime of validity of the simulations, we will need to 
bound ourselves to certain conditions. We need small string coupling (equation (36)) 

9s < 1 , (93) 



and small background fields (equation (37)) 

[m]^ 2/3 «l , [n]<?y 3 l. (94) 

To assure that the probe motion is slow enough, we need it to be much smaller in size than 
the length scales set by the background fields (equation (39)) 

[s} 2 [m] < 1 , [s] 2 [n] 2 < 1 . (95) 



Finally, to assure that the expansion of the DBI action is convergent, we need (equation (38)) 

- ,V3 « 1 



l s \9s 



[s} 2 [m]gt /3 « 1 . 



(96) 



The result of all these conditions is a region depicted in Figure [TT] We hence restrict again 
our simulations to the relevant parameter space. 



4.3 Case III: Falling into a black hole 

In this scenario, we would like to determine the time evolution of the probe DO branes as 
they approach the horizon of a four dimensional Schwarzschild black hole. For simplicity, we 
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Figure 11: The region of the parameter space for which our simulations can be trusted for 
Case II. Three planes bound it: one arising from convergence of the DBI expansion, another 
from requiring small background curvatures, and another from requiring a small probe size 
with respect to the background length scale to assure quasi-static evolution in the center of 
mass frame of the probe. All other bounds are weaker that the ones depicted in the input 
parameter range we consider. This is a semi-infinite region, with no bounds on the smaller 
ends of all three axes. 
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turn off background ffuctuations. In the center of mass frame of the probe, the metric looks 
like |al] 

^ = -(l + ^(, 2 + , 2 -2x 2 ))^ + ... (97) 

in Fermi normal coordinates. The x direction is the in-fall direction, while y and z are 
transverse. The remaining six directions of space are compactified - that is the probe's 
evolution is frozen in these directions, r is the location of the center of mass of the probe 
away from the black hole center (i.e. it is the Schwarzschild radial coordinate), and r is 
the time as measured by the in-falling probe. Hence, the dynamics is 3 + 1 dimensional. 
As the probe approaches the horizon, the metric evolves slowly compared to time scales 
associated with inter-probe dynamics. This is known as the adiabatic regime of the Fermi 
normal coordinates. The time dependence of r is given by 

r = ^(l + cosw) (98) 

where r is the initial radial position of the center of mass, and u is related to the local Fermi 
time t via 

Or 



T ^¥fe (w + sm ^ (100) 

We assume the in-fall starts with zero velocity. We then identify the background field 
parameters of our DBI action as 

-^ G tt,xx= ^- = M , --G tt , yy = --G tt , zz = -—^-=M 1 = ^, (101) 



where we have divided the three dimensional subspace into two subspaces as defined in ( 43 ) . 
The temperature of the black hole is given by 

with rh being the location of the black hole horizon. In dimensionless variables, we then 
have 

mi = _^ = 4 ^)V (103 ) 
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where t is the black hole temperature T in units of l/£. To find r as a function of time r, 
we start from (98) and (100) 



3/2 



87TtT 

with 

Expanding in small t/tq, we get 

r 



r 
r 



r 
r 



+ arccos 



2r 
r 



r < r . 



— ~ I 1 — Qirt t 
r 



2/3 



The time of flight, before the probe reaches the central singularity, is given by 



3/2 



1 



67T t 



We can write the background fields more conveniently as 

m x (0) 



mi(r) 



1 - 3 v /m 1 (0)r 



with the time of flight as 



T fi 



(104) 



(105) 



(106) 



(107) 



(108) 



(109) 



3vW0) 

where m!(0) is the initial tidal force felt by the probe at time equal to zero, away from the 
horizon at r = r$. And a similar expression with m\ — > —mo/2. The time to reach the 
horizon is then 



Tho 



1 



2nt 



(110) 



As input parameter, we get to specify the string coupling g s , the initial size of the probe 
(through s), mi(0), the back hole mass or equivalently horizon r^, and the initial position 
of the probe tq. For the analysis to be reliable, we need to satisfy the following conditions. 

We need small string coupling (36) 



g s < 1 ■ 
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And weak background fields (37) 



[m]g 2 / 3 « 1 



(112) 



This condition needs to be satisfied at all times. This means that, given a choice of valid 
rai(O), the simulation is still to break down as the probe approaches the black hole singularity. 
Hence, we need to keep track of this condition as a function of time to know when to stop 



the simulation. A proper DBI approximation requires (38) 

[*k 2/3 « i • 



(113) 



This latter condition is now trickier to handle: in cases I and II, the probe invariably collapsed 
to smaller sizes; hence, assuring that s is small enough at time zero was enough. In this 
case, we may see s grow in size as a function of time, as the probe approaches the black hole. 
Hence, we need to check that this condition is satisfied throughout time. Finally, these Fermi 
normal coordinates are reliable in the adiabatic regime. This requires the condition 

e = < 1 . 



(114) 



R is the size of the probe as defined in Once the probe becomes large enough to violate 
this condition, we need to stop the simulation. 

Our strategy goes as follows. Having experimented with the previous two cases, we know 
that we need a time of flight of the order of Tfu g ht ~ 500 with 5r = 0.01 to maintain control 
over the numerics of the simulations at the level of 10% error. This fixes mi(0) ~ 10 -6 . 
Hence, choosing any g s <C 1 is enough to handle the first two conditions (111) and (112). 



The third condition (113) is inconsequential since we can always adjust g s to be smaller 



(which makes the Planck length smaller) for given s. The fourth condition (114) is however 
one we need to keep a close eye on. It says that as the probe approaches the horizon, we 
need to make sure that the probe size R is much less than the horizon size. This may be a 
problem, since one school of thought is that the probe may spread to the size of the black 
hole as it crosses the horizon... If this is the case, the adiabatic approximation of the Fermi 
coordinates breaks down and we may need to revise the setup. 



5 Numerical techniques 

The simulations were implemented as a first order coupled Hamiltonian system using a 
fourth-order Runge-Kutta algorithm with a typical time increment of 5r = 0.01 - except 
when the thermal fluctuation timescale is smaller t > 100, where we reduce the size of the 
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time-step accordingly St ~ 1/t. The equations of motion are 

•i -A 



x" = p 1 
y a = q° 



f = [x j , [x\ x j ] ] + [y a , [x\ y a ] ] - 2m 1 x i - 3ine ijk [x j , x k ] 

q a = [y\[y\y b ]\ + [x\[y\x l ]\-2m 2 y a (115) 

with the constraint 

[x\f} + [y a ,q a } = 0. (116) 

The initial conditions are set up with all the p"s, y a 's, and g a 's equal to zero, while the x"s 
satisfying the SU(2) algebra 

x^st 1 , [a\ a 3 ] = ie ij k a k (117) 

in an iV x iV representation. In addition, small fluctuations were added to the x % matrix 
entries at the level of about 10% — 50%. Without these fluctuations, the evolution is found 



to be deterministic and oscillatory, as shown in Figure 12, Indeed, we use this to test the 
stability of our algorithm. We checked for both energy conservation and the constraint 
violation. This happens invariably as numerical errors accumulate over time, or when the 
probe starts exploding in size introducing large numbers in the simulation code that are 
more susceptible to computational errors. The simulation code was first implemented in 
Mathematica; and then, for efficiency purposes, recoded independently in objective-C. The 
two independent codes running on several different platforms, from Linux to Mac, were tested 
against each other and verified to yield the same results. Then the computation-intensive 
work was delegated to the objective-C version with the addition of thermal fluctuations. 
Single-precision floating numbers were used using the BLAS and LAPACK open source 
linear algebra libraries that take advantage of parallel-processing with multicore CPUs. A 
typical simulation for AT = 15 takes about five minutes to complete on 2 — 4 cores. At 
larger values of N, we ran into memory limitations and hence an associated disk-writing 
time penalty. Eventually, the code will be ported to GPU systems and is expected to run 
up to 100 times faster. For now, we restricted our simulations to AT < 100. 



With the fluctuations added, the constraint (116) is satisfied at the initial time. The 
system's evolution then assures that the constraint remains satisfied. We use the violation 
of the constraint to track numerical error and stop the simulation when the errors reach the 
10% level. 

The results shown in this work involved over 500 simulations with varied parameters. 
Mathematica was used to analyze the results and generate the plots. All code is available 
upon request. We also have developed a visual interface front-end to the code that allows 
realtime exploration of DO brane dynamics. 
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Figure 12: The time evolution of the probe when no initial fluctuations are added. We see 
that energy conservation and the constraint get violated by accumulated numerical errors 
after 40 — 50 natural oscillations of the system. For the constraint, we plot the average of 
the norms of the matrix entries in the constraint equation. 
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6 Conclusion and outlook 



We have found that DO brane probes have an internal mechanism, at strong coupling, for 
thermalizing and scrambling their initial conditions. The quartic term in the Hamiltonian 
is at the heart of this process. This was already determined for the BMN case in [2T]. Our 
results extend this conclusion to a large class of backgrounds. 

This process of thermalization is a quick one: the configuration reaches an equilibrium 
size and temperature after a few natural oscillations. The final temperature and size, as well 
as the time scale of thermalization, are independent of thermal background field fluctuations 
- a surprising outcome given that in most traditional statistical mechanical systems this 
process is intimately tied with the dynamics of the background thermal bath. Furthermore, 
background metric geometry on its own is not obviously correlated with this thermalization 
mechanism if the probe size is much smaller that the curvature length scale of the background; 
but some background RR fluxes do play a critical role. These observations suggest that, if 
a black hole geometry is to imprint its attributes onto the thermalization process of an in- 
falling probe, it can do so in one of two ways: (1) Have hair, perhaps along the line of the 
fuzzball proposal [H]. In our language, this corresponds to the effect of D2 brane flux on 
the thermal attributes of the probe; or (2) By expanding the in-falling probe to the order of 
the size of the horizon, a phenomenon we see a signature of in this work on the gauge theory 
side, and in previous work [25J on the gravity side. It is most likely that both mechanisms 
are concurrently at work in higher dimensions. 

Looking more specifically at the many scaling relations we extracted from the simula- 
tions, we found that the equilibrium size of the thermalized DO brane probe scales linearly 
with the number of DO branes N: the more DO branes, the bigger the equilibrium configu- 
ration, suggesting an incompressibility bound that these Planckian building blocks seem to 
be saturating. Furthermore, the larger this equilibrium size, the quicker the thermalization. 
And the equilibrium temperature of the probe increases with its equilibrium size, suggesting 
that perhaps we are not seeing the probe collapsing into an individual black hole - which 
would have an inverse relation between its temperature and its size. Instead, the probe may 
be merging with a large background black hole. We see a correlation between the size of 
the probe and background D2 flux akin to the dielectric effect of [20], but now in this more 
general dynamical scenario. The Matrix black hole picture of [221 |23j EI] has reinforced the 
idea that the black hole must have some order or structure, perhaps in the form of a fuzzy 
D-brane stretched at the horizon. While the effect of the D2 flux on the thermal properties 
of the probe such as its temperature is in tune with such a picture, we fail to see a signature 
of this order emerging within the DO brane probe itself. Perhaps this requires a delicate 
phase transition within the probe [25J which can be captured when more background fluxes 
are turn on; in particular, DO brane flux, which we did not consider, would give a magnetic 
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velocity dependent force that may play a crucial role. 

Intuitively, we can summarize the thermalization phenomenon we observe throughout our 
simulations as follows. At strong coupling, the DO brane probe collapses in size but maintains 
a coherent, dense, seemingly incompressible configuration of finite smaller size. And when 
the system is maintained in such a confined space long enough - a few natural oscillations of 
the probe - the non-linear quartic term in the Hamiltonian scrambles the information in the 
probe into a thermal distribution. The process is characterized by internal attributes of the 
probe, or external attributes of the background that effect the equilibrium size of the probe. 
This latter mechanism can arise from either background fluxes, or from the probe becoming 
large enough to experience the IR cutoff set by the background metric curvature. 

Through the works of [5] , it was proposed that black holes are the fastest scramblers 
in the universe. The associated scrambling time r t h is conjectured to scale as 

r th t p ~ log N, (118) 

where we have written the expression in terms of our notation (t p being the probe tem- 
perature). This is to be contrasted with normal matter scramblers associated with a time 
scale 

r th t p ~N a (119) 

where a is a dimension dependent number. We will now argue that our results are indeed 
consistent with the conjectures of [H |8]. To see this, we start from ((9j) and (10) and get 



Tthtp ~ s 3 N 2 . (120) 
s tunes the energy of the system, and hence is related to the dimensionless coupling defined 



in (27). Given we are probing a large N regime, we expect that the effective coupling would 



instead scale as 

N 



,2 

-3 



fyf = 13 ( 121 ; 



instead of e . e is the energy of the probe, and, at equilibrium, scales as e ~ s N 1 as argued 



in (33)-(35). 7 is close to 2 for a structure that is truly random. However, if there is SU{2) 
order at equilibrium - that is a membrane like structure to the DO branes - it would be close 
to 3 as seen earlier for the perfectly spherical initial configuration. Looking back at ((81), we 



see that, at equilibrium, we have 

Tr (aW 

~N 

which indeed implies that the equilibrium configuration is closer to one forming D2 branes 



R\ q = v ;/ J ~ s 2 N 2 (122) 



with 7 = 3. Putting things together with 7 = 3, we can express s in terms of g 2 ^ = N/e 3 
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in ( pOj ) 

Ttht P ~ — = (123) 

with no iV dependence left! The stronger the effective coupling, the shorter the thermaliza- 
tion or scrambling time, which makes sense. And we do not have a power law dependence on 
the degrees of freedom N. The conjecture of [HE] suggests a logiV dependence. We believe 
that to capture this leading behavior in N, we need a full quantum treatment of the problem. 
That all powers of iV cancel in this expression at the classical level is consistent with the 
conjecture. This also syncs well with our observation that D2 brane background flux affects 
the size of the configuration: if the background fields are to arise from other Matrix black 
hole DO branes interacting with the probe, a D2 brane order in the background Matrix black 
hole would be seeding this flux that is needed to correlate thermalization properties of the 
probe with those of the background black hole. 

For all of our simulations, a numerical instability developed after enough time had elapsed. 
It cannot be determined conclusively and with certainty whether this phenomenon is a real 
physical one or a pathology of the numerical evolution of the Hamiltonian system. For zero 
background fields, the system does have a flat direction that drives the system to larger 
entropically favored configurations - which in turn increases numerical instability. However, 
we see this numerical instability even when there is no flat direction in the system. While 
we determined scaling relations for the lifetime of the probe, we do not hence know whether 
to accord any physical significance to these. 

Our technique did eventually run against computational limitations and the truly large 
N behavior of the dynamics could not be extensively developed. Fortunately, this can easily 
be remedied given that our simulation code is parallelized and can be adopted relatively 
easily onto high performance GPU architecture. This would take us from 4 cores to around 
400 cores; that is a factor of 100 decrease in computation time. We plan to implement these 
technical improvements in the near future. 

Another direction for further exploration involves determining the effects of quantum 
back-reaction from the background fields onto DO brane dynamics. It has been suggested 
in the fuzzball proposal context that the many fuzzball geometries are in a quantum super- 
position [32]; hence a probe would be sensitive to quantum fluctuations of the background 
geometry instead of the thermal ones investigated in this work. This topic may be accessi- 
ble analytically using adiabacity techniques developed in the context of the study of Berry 
phases. We hope to report on this in the near future. 
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